#!/usr/bin/python
import sys
import os
from operator import itemgetter, attrgetter


# parameters.
pileup_dir = sys.argv[1]
agp_file = sys.argv[2]
out_file = sys.argv[3]

# build index of files.
files = {}
for f in os.listdir(pileup_dir):
	tmp = f.replace(".pileup","")
	files[tmp] = True
	

# determine pair with best coverage.
fin = open(agp_file,"rb")
lines = fin.readlines()
fin.close()


# loop over agp.
pcov = False
prev = ""
prev_scaf = ""
pair_cov = []
cnt = 0
for line in lines:
	# skip blanks.
	if line[0] == "#": continue
	tmp = line.split("\t")
	
	# get gap.
	if tmp[4] != "W":
		gapsz = int(tmp[5])
		continue
	
	# get contig info.
	ctg = tmp[5]
	scaf = tmp[0]
	
	# skip first.
	if prev == "":
		prev = ctg
		prev_scaf = scaf
		continue
		
	
	# ensure same scaffold.
	if prev_scaf == scaf:
	
		# read in profiles.
		curfs = ["%s/%s.pileup" % (pileup_dir, prev), "%s/%s.pileup" % (pileup_dir, ctg)]
		profs = []
		missing = False
		for i in range(2):
			# skip if already done.
			if i == 0 and pcov != False:
				profs.append(pcov)
				continue
				
			# check if pileup is available.
			if ctg not in files or prev not in files:
				missing = True
				break
				
			fin = open(curfs[i], "rb")
			lines = fin.readlines()
			fin.close()
			
			# calulate average coverage.
			tot = 0.0
			for line in lines:
				if line[0] == "#": continue
				
				# get coverage.
				tmp = line.split()
				for t in tmp:
					tot += float(t)
				
			# divide by size.
			cov = tot / float(len(lines))

			# record.
			profs.append(cov)		
			pcov = cov
		
		# ensure we weren't missing one.
		if missing == False:
			# get pair coverage.
			cov1 = profs[0]
			cov2 = profs[1]
			avg = (cov1 + cov2) / 2.0
			pair_cov.append((avg, prev, ctg, gapsz))
		
	# update contigs.
	prev = ctg
	prev_scaf = scaf
	cnt += 1
	
	if cnt > 1000: break
	print cnt

# print sorted.
tmp = sorted(pair_cov, key=itemgetter(0))
tmp.reverse()

fout = open(out_file,"wb")
for t in tmp:
	fout.write("%f\t%s\t%s\t%i\n" % t)
fout.close()

